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We investigate the thermodynamic properties of a field-induced supersolid phase in a 2D quantum 
antiferromagnet model. Using quantum Monte Carlo simulations, a very rich phase diagram is 
mapped out in the temperature - magnetic field plane, with an extended supersolid region where 
a diagonal (solid) order coexists with a finite XY spin stiffness (superfluid). The various quantum 
and thermal transitions out of the supersolid state are characterized. Experimental consequences in 
the context of field-induced magnetization plateau materials are briefly discussed. 
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Introduction — Since the early proposals of a possible su- 
persolid (SS) state of matter where both diagonal 
solid and off-diagonal superfluid (SF) long range orders 
coexist, the existence of this new exotic state of mat- 
ter at low temperature has been intensively discussed 
in 4 He [H, 0, 0, 0]- For lattice models, bosons have re- 
vealed an interesting tendency towards the formation of 
a stable (lattice) SS phase on frustrated networks like 
for instance the triangular lattice @, 0] • Perhaps a more 
promissing route for the search of such an exotic state of 
matter points towards gapped quantum magnets in an 
external magnetic field where, using a hard-core boson 
representation for triplet rungs a diagonal order is 
equivalent to a Bose solid which breaks the lattice trans- 
lational symmetry and an off-diagonal order represents a 
finite fraction of condensed bosons which breaks the U(l) 
symmetry. 

Several Mott insulators are known to exhibit field in- 
duced Bose-Einstein condensation of the triplets |T(| , 
as well as magnetization plateaus where an incompress- 
ible commensurate crystal of triplets is formed [ll|, Il2l . 
[f3j]. However, a direct continuous (second order) tran- 
sition between the U(l) symmetry breaking condensate 
and the charge density ordered plateau state would vio- 
late the Ginsburg-Landau paradigm. Therefore, one ex- 
pects either a first order transition, or an exotic decon- 
fincd critical point [l4| , or an intermediate region where 
both competing orders might coexist, namely a SS state. 
Thus, it is natural to consider quantum magnets as seri- 
ous candidates to achieve a SS state in the vicinity of a 
magnetization plateau. And indeed a SS phase has been 
recently reported for a spin-i bilayer system fTHj ] . 

So far, most of the theoretical investigations have fo- 
cused on the existence and the stability of the SS at 
T = 0, leaving unanswered the question of the nature of 
the quantum and thermal transitions out of the SS. Note 
however that the thermal melting of the SS via two dis- 
tinct transitions (Kozterlitz-Thouless and 3-state Potts) 
on the triangular lattice was reported in Ref. 0- 

In this Letter, using the stochastic series expansion 
(SSE) quantum Monte Carlo algorithm we explore 
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FIG. 1: (color online) Magnetic field - temperature phase di- 
agram of the bilayer system Eq. (flj) depicted on the top. The 
supersolid region (dark grey) is displayed together with its 
neighboring phases: superfluid, solid, paramagnet. The crit- 
ical lines belong to different universality classes, as indicated 
on the plot. The transition points (various symbols) result 
from quantum Monte Carlo simulations and finite size scaling 
analysis (see text and other plots). The thin dashed lines are 
constant B or T scans discussed below. 



the properties of the SS state of the quantum magnet 
model introduced in Ref. with emphasis on a precise 
characterization of the quantum and thermal phase tran- 
sitions out of the SS state. The resulting magnetic field 
- temperature phase diagram is shown in Fig. [TJ As we 
shall see below, both T = quantum phase transitions 
(SF-SS and SS-solid) as well as the thermal transition 
(SS-paramagnet) display special features of direct exper- 
imental relevance. In particular, the melting of the SS 
is a two-step process in contrast to standard transitions 
usually observed in quantum antifcrromagncts. 



Zero temperature properties — We consider the following 
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bilaycr spin-i XXZ Hamiltonian 

+j ± J2 s ^- s ^-wbB st ta , (i) 



where we set gps = 1 and choose the following pa- 
rameters J = 1, A = 3.3 and J± = 3.45, for which 
an extended zero temperature SS state was found in 
Ref. [l5|. Before discussing the finite temperature melt- 
ing of the SS, we first characterize the ground state (GS) 
properties, using large scale SSE simulations on finite 
systems (L x L x 2) at low enough temperature: we 
have checked that using T = (2L)" 1 insures that we 
focus on the GS properties. A diagonal order at a vec- 
tor q = {q x ,<ly,qz) will be signaled by the saturation 
to a finite value of the static structure factor per site 
5(q) = {(J2r ex P(*Q ' )/-^ 2 > where the sum over r 
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1,2 and 



runs over r x 

N = 2L 2 is the total number of sites on the bilaycr. For 
the superfluidity, in contrast to Ref. [l5| where the con- 
densate density was calculated, here we concentrate on 
the spin stiffness p s { of the U(l) off-diagonal ordering 
in the xy plane since this is the order parameter at fi- 
nite temperature. It is measured via the fluctuations of 
the winding number 17J] in the SSE simulation. Using 
a hard-core boson mapping for the triplets [1], a finite 
value for 5(q) corresponds to a commensurate charge 
density wave state (incompressible) for the bosons, and 
a non zero spin stiffness p s f must be regarded as a finite 
SF density, corresponding to a gapless (compressible) SF 
state. The total density of effective bosons is measured 
via the average magnetization per site m z = S z ) /N, 
shown in Fig. [2] (a) after an infinite size extrapolation 
with L = 8, 12, 16, 24, 32. Upon increasing the field, m z 
grows coutinuously up to a half saturated magnetization 
plateau at m z = 0.25. This gapped plateau state, called 
solid (tt, tt, 0), displays diagonal order at q = (tt, tt, 0), as 
shown in Fig. [5] (c). When decreasing the field from the 
plateau, there is a gapless SF phase with a finite stiffness 
p s f (see Fig. |5](b)). Interestingly the solid (tt, tt, 0) is also 
present in this SF regime (see Fig. H](c)), thus achieving 
a SS phase over a broad regime of magnetic field with 
both p s y£ and S(ir,ir,0) ^ 0. At lower fields, the solid 
is not stable anymore and a standard SF state with only 
p s ^ is recovered. 

The transition between the SF and the SS can already 
be observed by looking at the field-induced magnetization 
curve in Fig. [5] (a). Indeed, it displays a clear anomaly 
at the onset of the transition, as evidenced by a feature 
in the spin susceptibility xQ>) = dm z /dB (pointed by an 
arrow in the inset of Fig. [2] (a)). This quantum phase 
transition SF-SS, driven by the external field B, occurs 
for B ss = 7.175(5). Its universality class seems to be con- 
sistent with 3D (2+1) Ising as the critical decay of the 
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FIG. 2: Ground state properties of the model |T} focused on 
the superfluid, supersolid (light grey region) and solid regimes 
[scan 2]. The three main panels show infinite size extrapola- 
tions of quantum Monte Carlo estimates versus the external 
field B, obtained for T = (2L) _1 with 8 < L < 32. (a) The 
uniform magnetization m z grows up to a plateau m z = 0.25. 
Inset (a): spin susceptibility x(^) = dm z /dB. (b) The spin 
stiffness p s f decays in the supersolid regime and vanishes lin- 
early at B B = 7.52 like p sf ~ 0.757(B - B s ) (straight line). 
Upper inset (b): Universal crossing of p s L z at the critical field 
B s with 2 = 2. Lower inset (b): Data collapse of p s L z versus 
(B - B s )L} /v obtained with z = 2 and v = 1/2. (c) The 
solid order parameter S(ty, tt, 0) 1 ' 2 grows up for B > 7.175 
and saturate in the solid phase. Inset (c) Finite size scaling 
of iS(7r,7r, 0) plotted versus 1/L for various fields B indicated 
on the plot; the lines are quadratic fits in L _1 . 



structure factor S^tt, tt, 0) ~ L~ 2fi l v at B ss (shown as a 
straight line vs L^ 1 in the inset of Fig. [5] (c)) is com- 
patible with the 3D Ising value 2/3 /v = 1.037 The 
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transition between the SS and the solid occurs when the 
stiffness vanishes, at B s = 7.52(1), as shown in Fig.[2](b). 
We find an excellent agreement with the SF-insulator uni- 
versality class [Tjj for this transition since the finite size 
scaling analysis is fully combatible with the mean-field 
exponents z = 2,v = 1/2,(3 = 1/2, as demonstrated in 
the inset of Fig. 0(b). 

Finite temperature properties — We now turn to the finite 
temperature properties of the SS and discuss in more 
details the phase diagram presented in Fig. [TJ 
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FIG. 3: (color online) Top panel (a): KT transition SF- 
paramagnet at B — 7 [scan 3(a)]. Spin stiffness p s { versus T 
for various system sizes indicated by different symbols. The 
intersections with the dashed line at 2/ttT* are plotted in 
the inset versus L and fitted to the logarithmic finite size 
scaling form (see text) with Tkt = 0.295 and Lq = 0.346 
(dashed blue curve). Bottom panel (b): Ising transition solid- 
paramagnet at B = 7.9 [scan 3(b)]. Structure factor S(n, n, 0) 
versus T for system sizes 16 < L < 48. The critical temper- 
ature T c ~ 0.435 is obtained by getting a data collapse for 
s(7r,7r,0)L 2 ' 3/ ' y as a function of (T - T c )L 1/v with ,3 = 1/8 
and v — 1, as shown in the inset. 



At T = the SS phase is expected to display a true 
long range order for both diagonal and off-diagonal com- 
ponents, breaking respectively the lattice translation and 
the U(l) symmetries. On the other hand this is not true 



for T > since the continuous U(l) symmetry cannot 
be broken in 2D and only a quasi-long-range order is ex- 
pected in the xy plane. Nevertheless, the spin stiffness p s f 
can remain finite up to a KT temperature Tkt and the 
diagonal order can persist at T > since it only breaks 
a discrete symmetry. Therefore, a persistence of the SS 
state at finite temperature is expected, but whether SF 
and solid components disappear simultaneously or at dif- 
ferent temperatures is not straightforward [Icj . As shown 
in Fig. [TJ we actually find that upon heating the melting 
of the SS occurs in two distinct steps: the SF density is 
destroyed via a KT transition while the solid order melts 
via an Ising transition at T c ^ Tkt- 

These critical properties of the SS are in fact intimately 
tied to the critical behavior of the surrounding SF and 
solid orders which are decoupled, and remain decoupled 
also in the SS regime. Indeed, as shown in Fig. [TJ the 
KT transition line between the SF and the paramagnct 
continues in the SS to eventually vanish at the SF-solid 
(T = 0) quantum critical point. Similarly, the Ising tran- 
sition line between the solid and the paramagnet is also 
prolongated down to the SF-SS quantum critical point. 

The KT and 2D Ising natures of respectively the SF- 
paramagnet and solid-paramagnet transitions are illus- 
trated in Fig. [3] where the scans 3(a) and 3(b) at B = 7 
and B = 7.9 are shown. In Fig. [3] (a), the disappear- 
ance of the 2D SF density p s { at Tkt is displayed for 
L = 8, 16, 24, 32. The KT transition, identified by a uni- 
versal jump for an infinite system /? s f(2KT) = —Tkt, su f~ 
fers from logarithmic finite size corrections and the in- 
tersection occurs at T*(L) = T KT [1 + 1/(2 ln(L/L ))], as 
shown in the inset of Fig. [3] (a). 

Upon heating, the (ir, n, 0) solid melts via an Ising 
transition into a paramagnct at a critical temperature 
T c , as shown in Fig. [3] (b). Using the finite size scaling 
hypothesis for the squared order parameter S(tt, vt, 0) = 
L-VJ/ujr j-( T _ T c )L l > v ] , the Ising universality class is 

checked by getting a data collapse for S(tt, n, 0)L 2 P/ u as 
a function of (T - T C )L 1 ^ for L = 8, ... ,48, with the 
Ising exponents (3 = 1/8 and v = 1 (inset of Fig. [3] (b)). 

As already stated, the melting of the SS occurs in two 
steps, as illustrated in Fig. 2] for B = 7.3 (scan 4 of the 
phase diagram Fig. [T]) . The universality classes (KT for 
the SF and Ising for the solid) are checked and the critical 
temperatures Tkt and T c are estimated. We notice that 
the logarithmic effects at the KT transition are more pro- 
nounced than for the SF to paramagnet transition since 
we found a length scale Lq for the logarithmic correc- 
tion much bigger than for the SF-paramagnet transition. 
While upon heating the transitions SS — > Solid — > Para- 
magnet (scan 4 for instance) are clearly indentified as 
KT and Ising 2D, it is much harder to conclude on the 
nature of the finite T transition SS — » SF. Indeed, even if 
2D Ising is also expected, finite size effects turned out to 
be very large and it would require a huge computational 
effort to settle this issue. 
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FIG. 4: (color online) Two different finite temperature tran- 
sitions for the destruction of the SS phase at B = 7.3 [scan 4]. 
Top panel: (2D Ising transition with /3 = 1/8, v = 1) struc- 
ture factor plotted 5(71", n, versus T for system sizes 
16 < L < 48. Inset: Finite size scaling with T c = 0.375. Bot- 
tom panel: (KT transition) spin stiffness p s f plotted versus 
T for the same system sizes as in the top panel. The dashed 
line 2T/n intersects the data points for T*(L), plotted in the 
inset vs l/ln(L/Lo) and fit to the logarithmic form (see text) 
with T KT = 0.149 and L = 5.25 (dashed blue line). 



Conclusions — We have confirmed the existence of the 
analog of a supersolid phase in a spin-1/2 bilayer, and 
we have fully characterized its temperature behaviour, 
which generically consists of two well separated phase 
transitions in the Kosterlitz-Thouless and Ising univer- 
sality classes. We have also shown that a clear magneti- 
zation anomaly is expected upon entering the supersolid 
phase from the supcrfluid one, and that the stiffness van- 
ishes linearly at the edge of the plateau, in agreement 
with the superfluid-insulator universality class. These 
properties are expected to be relevant for other geome- 
tries, in particular for systems where plateaus are pro- 
duced by frustration rather then anisotropic coupling, as 
recently agued in the context of a spin-1 model [2l[. As 
such, they should be useful in the search for supersolid 
phases in the neighborhood of magnetization plateaus 
observed in some quan tum antiferromagnetic compounds 
like SrCu 2 (B0 3 ) 2 lH. 

Transposed for bosons, the present results open a new 
route to supcrsolidity. Indeed, to first order in 1/Jj_, the 
effective hard-core boson model is the simple t — V model 
on the square lattice, which is well known to have phase 
separation but no supersolid phase Q . So the supersolid 
phase has to emerge from higher terms, such as correlated 
hoppings, which were recently shown to induce a pairing 
phase at low density Q . 

Finally, from a more theoretical point of view, we have 



shown that the T = quantum phase transition from the 
supcrfluid to the supersolid is consistent with a 3D Ising 
universality class. However, to which extent our numer- 
ical data are consistent with the general expectation for 
a transition between a superfluid and a checkerboard 
supersolid [24[ deserves further investigation. 
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